New gene signature from the dominant infiltration immune cell type in osteosarcoma predicts overall survival

The immune microenvironment of osteosarcoma (OS) has been reported to play an important role in disease progression and prognosis. However, owing to tumor heterogeneity, it is not ideal to predict OS prognosis by examining only infiltrating immune cells. This work aimed to build a prognostic gene signature based on similarities in the immune microenvironments of OS patients. Public datasets were used to examine the correlated genes, and the most consistent dominant infiltrating immune cell type was identified. The LASSO Cox regression model was used to establish a multiple-gene risk prediction signature. A nine-gene prognostic signature was generated from the correlated genes for M0 macrophages and then proven to be effective and reliable in validation cohorts. Signature comparison indicated the priority of the signature. Multivariate Cox regression models indicated that the signature risk score is an independent prognostic factor for OS patients regardless of the Huvos grade in all datasets. In addition, the results of the association between the signature risk score and chemotherapy sensitivity also showed that there was no significant difference in the sensitivity of any drugs between the low- and high-risk groups. A GSEA of GO and KEGG pathways found that antigen processing- and presentation-related biological functions and olfactory transduction receptor signaling pathways have important roles in signature functioning. Our findings showed that M0 macrophages were the dominant infiltrating immune cell type in OS and that the new gene signature is a promising prognostic model for OS patients.


Data collection
Gene chip expression data and the clinical information of OS patients from the TARGET database were obtained from the website https:// ocg.cancer.gov/ progr ams/ target/ data-matrix on May 1, 2022.Osteosarcoma was used as the only key word to select suitable datasets published up to May 2022 for our study on the Gene Expression Omnibus (GEO) datasets website (https:// www.ncbi.nlm.nih.gov/ gds/).Then, the datasets were searched by the following inclusion criteria: (1) the diagnoses were pathologically confirmed; (2) the dataset had Homo sapiens samples; (3) the platform contained whole-genome information; (4) the dataset included patient survival data; and (5) the sample size of the datasets was more than 30.After that, only the GSE21257, GSE16091, and GSE39055 datasets were included in our study.Therefore, the normalized mRNA expression data and survival data of the above datasets were downloaded from GEO (https:// www.ncbi.nlm.nih.gov/ geo/ query/ acc.cgi? acc= GSE21 257, https:// www.ncbi.nlm.nih.gov/ geo/ query/ acc.cgi? acc= GSE16 091, https:// www.ncbi.nlm.nih.gov/ geo/ query/ acc.cgi? acc= GSE39 055).
If the same patient provided two or more tumor samples to those datasets, only the data corresponding to the primary lesion were selected according to the sample numbers.All patients were included in our study.Gene symbols were extracted from the provided documents from those dataset websites.We conducted signal intensity normalization across the arrays of the above datasets using the normalizeBetweenArrays function from the limma package in R software.

Gene signature construction and validation
Gene signatures were generated by inputting M0 macrophage-correlated genes into the LASSO Cox regression model using the TARGET dataset.The "glmnet" package in R was used to complete the regression process.After the genes were selected, a multivariate Cox analysis was used to calculate the corresponding coefficients.The signature risk score was calculated by the sum of the products of each gene and its corresponding coefficients as follows: score = (CPE*0.7021)+ (FHL5*− 0.9794) + (GBP1*− 0.4718) + (GNLY*0.9393)+ (GPR82*− 1.8384) + (IL 18RAP*− 1.9630) + (LILRA2*0.4310)+ (NDRG4*− 0.7661) + (PLB1*− 0.6822).Moreover, the risk scores for the GSE21257 and GSE39055 datasets were also calculated using the same method with the coefficients derived from the TCGA dataset for validation.OS patients from each dataset were divided into two groups (low-and high-risk groups) according to the median value of the risk score to maximize the statistical power and provide the most reliable results based on a fixed sample size in those databases 33 .An overall survival analysis between the two groups was conducted using Kaplan-Meier curves with Wilcoxon log-rank tests.Multivariate Cox regression models were applied to test the independent predictive value of our signature.The prognostic accuracy of the risk score in different datasets was determined using the time-dependent receiver operating characteristic curve (ROC) and incident/dynamic (I/D) area under the curve (AUC) through the "timeROC" and "risksetROC" packages in R separately 34 .

Gene signature comparison
We screened the studies using the keyword "gene signature prognosis osteosarcoma" on PubMed.In addition, the inclusion criteria were as follows: (1) the journal impact factor was more than 5 and had a good reputation; (2) the online publication date range was from March 01, 2020, to March 01, 2022; and (3) the gene signature was constructed by messenger RNA to be easily validated in other datasets.After the studies were selected, the gene signatures were recalculated by a multivariate Cox proportional hazards model using data in a training cohort.Kaplan-Meier curves, time-dependent ROCs and I/D AUCs of the gene signatures using validation cohorts were chosen to compare the prognostic accuracy.

Gene set enrichment analysis (GSEA)
GSEA was launched to study the biological functions of our gene signature in the TARGET, GSE21257, and GSE39055 datasets using the "clusterProfiler" R package based on the GO and KEGG analyses individually.The GO functions and KEGG pathways with P values less than 0.05 in each dataset were merged.Finally, the biological functions with consistent positive or negative values of enrichmentScore were reserved to find the real difference of GO functions and KEGG pathways between low-and high-risk groups discriminated by our signatures.

Correlation between signature risk score and chemotherapy sensitivity
The IC50 values of chemotherapy drugs (bleomycin, cisplatin, doxorubicin, etoposide, and methotrexate) in each osteosarcoma sample in the TARGET, GSE21257, and GSE39055 datasets were accessed through the "pRRophetic" package in R 35 , which was built based on the Genomics of Drug Sensitivity in Cancer database (www.cance rRxge ne.org) 36 .Then, the samples were classified into low-or high-risk groups by the median of our gene signature risk score.Finally, the IC50 values of chemotherapy drugs between different groups were analyzed by the Mann-Whitney rank test.A P value < 0.05 was defined as statistically significant.

Statistical analysis
We conducted all statistical analyses with R version 4.0.5 (R Foundation for Statistical Computing, Vienna, Austria) and GraphPad Prism 6.01 (GraphPad Software, Inc., San Diego, CA, USA).OS patients were divided into two groups by the median of our signature risk score in the TARGET, GSE21257, and GSE39055 datasets.Violin plots were performed in Hiplot (https:// hiplot.com.cn), a comprehensive web platform for scientific data visualization.The comparison of 22 immune cells, immunotherapy-targeted genes, ESTIMATE scores, stromal scores, and immune scores between low-and high-risk groups were analyzed using a Mann-Whitney rank test.Unless otherwise specified, a P value < 0.05 was defined as statistically significant.

Ethics approval
Our study is based on open-source data (TRGET and GEO).Ethical review and approval were not required for the study on human participants in accordance with the local legislation and institutional requirements.www.nature.com/scientificreports/

Consent to participate
All methods were carried out in accordance with the relevant guidelines and regulations.

Clinical characteristics of the studied datasets
Figure 1 shows the flowchart of this study.The clinical characteristics of OS patients in the TARGET, GSE21257, GSE16091, and GSE39055 datasets are shown in Table 1.In the TARGET and GSE21257 datasets, the age, sex, tumor location, and metastatic status at diagnosis were similar.The GSE21257 dataset contained the most detailed clinical information except for radiotherapy treatment and had a similar Huvos grade distribution as the GSE39055 dataset.GSEGSE16091 only included survival information.The median overall time and survival curves (Supplemental Figure 1) were close to each other in all datasets.

Estimation of tumor immune cell infiltration in osteosarcoma
To obtain the landscape of tumor immune cell infiltration, we applied the CIBERSORT algorithm to the selected datasets.As shown in Fig. 2, the major immune cells in the TARGET and GSE16091 datasets were CD8 + T cells and M0 and M2 macrophages, while the major immune cells in the GSE21257 dataset were M0 and M2 macrophages.For the GSE39055 dataset, M0 macrophages were the most dominant infiltrating immune cells.Other immune cell infiltration levels were quite low in all datasets.In general, M0 macrophages were the dominant immune cells across all datasets.As shown in Supplemental Figure 1, the survival difference between each dataset was not significant, and we believe that M0 macrophages are influential on the overall survival of osteosarcoma.Therefore, the M0 macrophage-associated genes were screened by the methods described in the "Methods and Materials", and the selected genes are shown in Supplemental Table 1.

Construction and validation of prognostic gene signature for osteosarcoma
The above selected genes were used to build a prognostic signature in OS patients through LASSO Cox regression analysis by the TARGET dataset with 86 patients as the discovery cohort (Fig. 3A).An optimal 9-gene prognostic signature was made (Fig. 3B,C)).The biological functions and the coefficients of signature genes are shown in Table 2.The signature risk scores are equal to the sum of the product of the expression value and coefficient of each gene.We chose only the GSE21257 and GSE39055 datasets as the validation cohort to measure the prognostic value of the signature-based risk score, and the GSE16091 dataset was excluded owing to the absence of gene expression data for GPR82 and PLB1.We used the log-rank test to study the association of the risk score with survival in the TARGET dataset.As expected, high-risk patients had significantly shorter survival than low-risk patients (log-rank test, P < 0.0001) (Fig. 3D).The same tendency was confirmed in the GSE21257 and GSE39055 www.nature.com/scientificreports/datasets (log-rank test, P = 0.034 and P = 0.031) (Fig. 4A,B).In addition, the AUCs of the time-dependent ROCs including the signature risk score and selected clinical factors (age, sex, and Huvos grade) reached 0.793, 0.885, and 0.839 for 1-, 3-, and 5-year OS in the TARGET dataset, respectively, which were higher than the AUCs of the time-dependent ROCs including only clinical factors (Fig. 5A, B, C).Similar results were also observed in the GSE21257 and GSE39055 datasets (Fig. 5E,F,G,I,J,K).We compared the C index of I/D AUC to assess the discriminative accuracy of the prediction model with and without the signature risk score.The results showed that the C index increased to 0.831 (95% CI, 0.750-0.868P < 0.0001), 0.690 (95% CI, 0.527-0.754;P < 0.0001), and 0.775 (95% CI, 0.584-0.840;P < 0.0001) when the signature risk score was added to the prognostic model in the TARGET, GSE21257, and GSE39055 datasets (Fig. 5D,H,L).

Investigation into the differences in immune cells, immunotherapy-targeted genes, and immune scores between high-and low-risk groups
As the signature risk score reflects the tumor immune activity of OS, we further investigated the differences in immune cells, immunotherapy-targeted genes, and immune scores between different risk groups from all datasets.The immune cell infiltration results were inconsistent between the three datasets (Supplemental Figure 2).M0 macrophages were found to be higher in high-risk groups than in low-risk groups in the TARGET and GSE21257 datasets (Fig. 6A, B).Activated dendritic cells were found to be lower in the high-risk groups than in the low-risk groups in the GSE21257 and GSE35099 datasets (Fig. 6B, C).For the results of six hot immunotherapy-targeted genes and immune scores, no significant differences were found in the GSE35099 datasets (6F, I), which may be due to its relatively small sample size.TIM-3, LAG-3 and all immune scores were all found to be lower in the high-risk groups than in the low-risk groups in the TARGET and GSE21257 datasets (Fig. 6D,E,G,H).In summary, our signature risk score may reflect the immune activity of osteosarcoma.In addition, LAG3 may be the most effective immunotherapy target because its expression was relatively high in all datasets (Fig. 6D,E,F).Owing to the low incidence rate of OS, the OS datasets did not contain large sample sizes.Therefore, adjustments for multiple testing were not made to reserve important clues.

Comparison of gene signatures
Thirteen previously reported studies were included in the comparison after the selection process, which are shown in Table 3. Kaplan-Meier curves (Fig. 7) and univariate Cox models (Table 4) showed that only the difference in overall survival for the nine-gene signature model was significant in the training and testing cohorts.
Additionally, the AUCs of the time-dependent ROCs of the signature risk score were also found to be greater  than the AUCs of the time-dependent ROCs of previous signatures in the TARGET and GSE39055 datasets (Fig. 8).The results showed that our nine-gene prognostic signature is more robust than previously reported gene signatures.

Identification of signature-associated biological functions and pathways in OS
The tumor sample tissues of the three datasets were dichotomized into high-and low-risk groups according to the median of the nine-gene signature risk score.A GSEA was performed to identify the signature-associated biological functions and signaling pathways, and then the results from all datasets were merged to attain reliable results.We found that the olfactory transduction receptor signaling pathway and its related biological functions were downregulated in high-risk groups compared to low-risk groups (Fig. 9).Antigen processing-and presentation-related biological functions were also found to be upregulated in high-risk groups compared with low-risk groups.However, none of the related signaling pathways were found to be significantly correlated with different risk groups.These results indicate that the olfactory transduction receptor signaling pathway may play a role in the difference between signature-predicted outcomes.

Signature risk scores and chemotherapy sensitivity
To investigate whether the immune-correlated gene signature is an independent prognostic factor of OS patients, we need to adjust it based on previously reported important prognostic factors, including histologic response to chemotherapy and presence of metastases.Supplemental Figure 3 shows the results and indicates that our signature is an independent OS prognostic factor.However, not all target patients received neoadjuvant chemotherapy; thus, this group of patients lacked Huvos grade information.To validate the independence of our signature, we explored the relationship between the risk score and chemotherapy sensitivity.IC50 was calculated using the "pRRophetic" R package to predict the treatment response to chemotherapy drugs (bleomycin, cisplatin, doxorubicin, etoposide, and methotrexate).As we anticipated, there was no significant difference in the sensitivity of all drugs in the low-and high-risk groups (Fig. 10).The results are consistent with the multivariate Cox analysis.GSEA also did not reveal any biological functions or pathways responsible for chemotherapy resistance.Therefore, the overall survival difference predicted by our signature is probably due to the various backgrounds of the immune microenvironment of OS.

Discussion
Although OS is the most common primary bone cancer in children and young adults, it is a very rare cancer, with approximately 400 new cases diagnosed annually in the USA 37 .The incidence peaks in adolescence and in old age 38 .The most common early symptom is ostalgia, which is easily confused with growing pain.The most common location of OS is in the metaphysis around the knee joint, followed by the proximal tibia and humerus.www.nature.com/scientificreports/ The majority of patients are therefore diagnosed with localized disease.Among the patients in the progression stage, the lungs and other bones are the most common metastatic locations 39 .The prognostic factors include tumor site, tumor size, tumor resectability, histologic response to chemotherapy, and presence of metastases [40][41][42] .Chemotherapy consisting of high-dose methotrexate and other drugs before and after definitive surgery has been used as the standard treatment strategy for localized osteosarcoma over the past decades [43][44][45][46] .In patients with relapsed and/or metastatic osteosarcoma, surgery and chemotherapy have some effect 3,[47][48][49] .However, the prognosis of those patients was a median event-free survival of less than 4 months, which is unsatisfactory 50 .An increased understanding of the disease from well-annotated tissue banks through newly developed technologies has revealed its heterogeneity and molecular aberrations, which offer new insight into targeted therapy and immunotherapy.
Recently, immunotherapy has facilitated a revolution in many kinds of solid malignant tumor treatments.Although its efficiency is not so high, it can lead to considerable survival benefits among responders 51 .Osteosarcoma humanized in vitro and humanized mouse models with similar tumor microenvironments in several previous studies have indicated that immunotherapy targeting the tumor microenvironment may be a promising strategy for osteosarcoma treatment [11][12][13] .Studies aiming at finding the precise immune response population have also suggested that the tumor immune microenvironment is the critical factor 52 .Therefore, it is essential to better understand and characterize the immune status of OS for the precise prediction of the immune response and the development of immunotherapy.The OS immune microenvironment consists of a network of immune cells that function as ideal grounds for tumor proliferation and progression 15,53 .
Our study first found that M0 macrophages were the dominant immune cells in the TARGET and GEO databases using the CIBERSORT algorithm.However, the direct association between M0 macrophage infiltration and the overall survival of OS patients was not consistent in different independent datasets (Supplemental Figure 4).Based on the hypothesis that the dominant immune cells play the most important role in the OS microenvironment, we screened for genes associated with M0 macrophages and built a new nine-gene signature for OS prognosis using clinical information from the TARGET dataset.In addition, the prognostic signature was proven to be effective and reliable in validation cohorts consisting of two independent GEO datasets.Previously published prognostic gene signatures were compared with the newly discovered gene signature using publicly available datasets.The results showed that no signatures besides ours offer a precise prediction value of OS survival in all datasets, indicating the priority of the nine-gene signature.Multivariate Cox regression models indicate that the signature risk score is an independent prognostic factor for OS patients regardless of age, sex, or Huvos grade in all datasets.Furthermore, the association between the signature risk score and chemotherapy sensitivity also showed that there was no significant difference in the sensitivity of all drugs between the low-and high-risk groups.Therefore, the different outcomes predicted by our signature are probably due to the various immune microenvironments of OS.To our knowledge, this is the first time a study has been conducted on the dominant infiltrating immune cell profiles in several datasets to find similarities between OS tumor immune microenvironments.As we expected, a robust gene signature was generated based on the discovered similarities.
M0 macrophages are slightly elongated unstimulated macrophages, which are considered to be theoretically inactivated 54 .Previous studies have suggested that M0 macrophages are associated with poor prognosis and metastatic disease in lung adenocarcinoma 55,56 , hepatocellular carcinoma 57 , osteosarcoma 53 , pancreatic adenocarcinoma 58 , melanoma 59 , colorectal cancer 60 , gastric cancer 61 , and glioblastoma 62 .Recent studies have shown that new tumor-infiltrated M0 macrophages induce pancreatic cancer cell death via TNF-α secretion, while M1, M2, and tumor-associated macrophages do not harbor antitumorigenic activities 63 .However, mesenchymal stem cells could activate STAT6 and induce M2 polarization through disease progression, leading to anti-inflammatory functions 64 , which may explain the relationship between high M0 macrophage infiltration and poor prognosis.Our prognostic signature consisted of nine genes (Table 2): CPE, FHL5, GBP1, GNLY, GPR82, IL18RAP, LILRA2, NDRG4, and PLB1.In our research, CPE was determined to be an unfavorable biomarker for OS prognosis, while PLB1 exhibited protective effects in OS (Supplemental Figure 5).CPE encodes carboxypeptidase E, which belongs to the carboxypeptidase family and is reported to be involved in the biosynthesis of hormone and neurotransmitter peptides and to play nonenzymatic roles in the endocrine and nervous systems 65 .However, previous studies have also proven that CPE plays various roles in various cancers 66 .Its N-terminal truncated protein has been found to be expressed in multiple cancers, regulate metastatic gene expression, and encompass different signaling pathways 67 .CPE has been reported to reduce aerobic glycolysis and migration in glioblastoma cells 68 .In addition, CPE may promote migration and invasiveness in various other cancers, such as lung www.nature.com/scientificreports/cancer 69 , pancreatic cancer 70 , and cervical cancer 71 .In osteosarcoma, Fan et al. found that CPE-variant proteins increase invasiveness and promote epithelial-mesenchymal transition through the activation of the Wnt pathway in OS cells 72 .CPE also modulates immunity.Bar et al. found that it regulates immune homeostasis and inhibits inflammation 73 .Sanjay also reported that it enhances the innate immunity of the male reproductive tract 74 .Therefore, whether CPE plays a role in cancer immunity needs to be discussed further.PLB1 has seldom been studied in cancer.Its rearrangement with ALK has been found in lung cancer 75 , which indicates a sensitivity to targeted therapy 76 .Lin et al. reported that overexpressed PLB1 antigens cause high infiltration of immune cells and a favorable prognosis in glioblastoma patients 77 .PLB1 is commonly discussed in yeast and was found to be essential for the survival and virulence of Cryptococcus 78 and Candida albicans 79 .
Recently, with the appearance of advanced bioinformatics algorithms, numerous gene signatures have been built from public online datasets.However, the universalism of those signatures was largely undetected or only detected in one independent cohort.Regarding the limitations of previous findings, the nine-gene signatures were summarized from the generality of several datasets.To determine the priority of our hypothesis, we included previous potential survival-related OS signatures from recent studies published in journals with a high influence and good reputation for comparison [80][81][82][83][84][85][86][87][88][89][90][91][92] .As we expected, after validating them in the training and validation cohorts, the nine-gene signature was the only signature showing a consistent predictive value for OS survival.Among the previous studies, the signatures by Fan et al. and Feng et al. could accurately predict the prognosis of OS in the TARGET and GSE21257 datasets.However, in the GSE39055 dataset, the lines of different groups predicted by Feng's signature in the K-M plot were very close, while the survival difference of groups predicted by Fan's signature showed a reverse tendency compared with other datasets.
A GSEA of GO and KEGG pathways found that antigen processing-and presentation-related biological functions and olfactory transduction receptor signaling pathways play important roles in signature functioning.Numerous studies have proven that antigen processing-and presentation-related biological functions are key functions in the cancer immune response process and are mainly carried out by dendritic cells, macrophages and B cells [93][94][95] .These results indicate the important role of the three antigen-presenting cells and cancer immunity in the prognosis of OS, which is consistent with our hypothesis.Regarding the olfactory transduction receptor signaling pathway, a previous study reported that olfactory receptors provide innate and adaptive immune responses during the virus entry process 96 .Orecchioni et al. recently found that immune cells, such as vascular macrophages, express olfactory receptors, which induce interleukin-1β secretion, leading to inflammation 97 .Nevertheless, there have not been any studies that have reported any association between the olfactory transduction receptor signaling pathway and cancer immunity thus far.
Our analysis also showed the homogeneity and heterogeneity of the osteosarcoma immune microenvironment.M0 macrophages were the dominant infiltrating immune cell type in all datasets.In addition, the infiltration of M2 macrophages, which have been reported to play an immune repressive role in cancer 98 , was similar to that of M0 macrophages in the three datasets, except for the GSE39055 dataset.There were also several infiltrated CD8 + T cells in the TARGET and GSE16091 datasets.Other immune cell infiltration levels were quite low.The above results may partly explain the unsatisfactory response to current immunotherapy in osteosarcoma.Most of the immune microenvironment of osteosarcoma tissues consists of immune repressive M2 macrophages or nonfunctional M0 macrophages.Patients with tumor tissues infiltrated with several CD8 + T cells may receive survival benefits from immunotherapy.Given the large proportion of M0 and M2 macrophages in osteosarcoma immune cells, inducing M0 or M2 macrophages toward the M1 phenotype to promote the antitumor immune response may be a promising treatment strategy for OS patients.There are several limitations in the present study.First, the prognostic signature was constructed from public online retrospective data with relatively small sample sizes.It may be improved by future OS studies with larger sample sizes.Second, the predictive value of the signature needs to be confirmed by future prospective studies.Third, the mechanism and function of a high infiltration of M0 macrophages in the OS microenvironment is still unclear and needs to be clarified in the future.Fourth, many prognostic factors were unavailable in the datasets used, and as such, the independence of the signature could not be fully determined.In addition, there are no wet-lab experimental data supporting the roles of signature genes and the olfactory transduction receptor

Conclusions
In the present study, we built a new and robust nine-gene signature based on the hypothesis that the dominant infiltrated immune cells play the most important role in OS progression.The newly defined gene signature was found to be significantly associated with OS prognosis in all datasets.In addition, we proved the advantage of the signature by comparing it to previously published signatures.Antigen processing-and presentation-related biological functions and the olfactory transduction receptor signaling pathway were found to be associated with the signature risk score.The potential role and mechanism of the olfactory transduction receptor signaling pathway and M0 macrophages in OS should be evaluated in the future. https://doi.org/10.1038/s41598-023-45566-6

Figure 1 .
Figure 1.A schematic diagram of dataset analysis.
Figure3Eand Fig.4C,Dshow the distribution of signature gene expression profiles and survival statuses in different risk score groups for OS in all datasets.Figure3Fand Fig.4E,Fshow that the high-risk score was significantly associated with a risk of increased mortality in OS patients in the multivariate Cox regression models after adjusting for age, sex, and Huvos grade in all datasets (HR = 51.871,3.204 and 6.663, 95% CI = 6.898-391.183,1.066-9.637and 1.227-36.171,P < 0.001, = 0.038 and = 0.028, respectively), which indicates that the signature risk score is an independent prognostic factor for OS patients.

Figure 2 .
Figure 2. Estimation of tumor immune cell infiltration in osteosarcoma by the CIBERSORT algorithm.The percentage of infiltration of 22 immune cell types in the TARGET (A), GSE21257 (B), GSE16091 (C), and GSE39055 (D) datasets showed that M0 macrophages were the dominant immune cell type.

Figure 3 .
Figure 3. Construction of a prognostic gene signature.(A) Venn diagram of the M0 macrophage-associated genes in the TARGET and GSE21257 datasets.The numbers in each area represent the gene numbers in each group.(B) Cross-validation for tuning parameter screening upon LASSO regression analysis.(C) Screening of the optimal parameter (lambda) at which the vertical lines were drawn.(D) Kaplan-Meier overall survival analysis of the gene signature risk score in OS of the TARGET dataset.(E) Distribution of signature gene expression profiles along with survival status in different signature risk score groups in TARGET datasets.(F) Forest plot showing the results of multiple factors in the Cox regression analysis of the gene signature risk score with other clinical characteristics in OS of the TARGET dataset.

Figure 4 .
Figure 4. Validation of the nine-gene prognostic signature.Kaplan-Meier overall survival analysis of the gene signature risk score in OS from the GSE21257 (A) and GSE39055 (B) datasets.Distribution of signature gene expression profiles along with survival status in different signature risk score groups in the GSE21257 (C) and GSE39055 (D) datasets.Forest plot shows the results of multiple factors in the Cox regression analysis of the gene signature risk score with other clinical characteristics in overall survival from the GSE21257 (E) and GSE39055 (F) datasets.

Figure 5 .
Figure 5. Prognostic values of the signature risk score in the training and validation cohorts.Time-dependent ROCs at 1 (A), 3 (B), and 5 (C) years and I/D AUC (D) for OS in the TARGET dataset.Time-dependent ROCs at 1 (E), 3 (F), and 5 (G) years and I/D AUC (H) for OS in the GSE21257 dataset.Time-dependent ROCs at 1 (I), 3 (J), and 5 (K) years and I/D AUC (L) for OS in the GSE39055 dataset.ROC, receiver operating characteristic curve; I/D AUC, incident/dynamic area under the curve.

Figure 6 .
Figure 6.Association between immune cell infiltration levels, immunotherapy-targeted gene expression, ESTIMATE immune scores, and signature risk score of OS.The immune cell infiltration types found to be significantly associated with the signature risk score of OS in the TARGET (A), GSE21257 (B), and GSE39055 (C) datasets.Association between immunotherapy-targeted gene expression and the signature risk score of OS in the TARGET (D), GSE21257 (E), and GSE39055 (F) datasets.Association between ESTIMATE immune scores and signature risk scores of OS in the TARGET (G), GSE21257 (H), and GSE39055 (I) datasets.

Figure 7 .
Figure 7. Comparisons of the gene signature with previously published gene signatures in the TARGET, GSE21257, and GSE39055 datasets using the Kaplan-Meier estimator.The results are in bold and considered significant if the P value < 0.05.

Figure 8 .
Figure 8. Comparisons of the gene signature with previously published gene signatures in the TARGET, GSE21257, and GSE39055 datasets using time-dependent ROC curves.

Figure 9 .
Figure 9. Gene set enrichment analysis of GO and KEGG pathways in OS between different signature risk score groups in the training and validation cohorts.The results of GO functions between different signature risk score groups in the TARGET (A), GSE21257 (B), and GSE39055 (C) datasets.The results of KEGG pathway analysis between different signature risk score groups in the TARGET (D), GSE21257 (E), and GSE39055 (F) datasets.

Figure 10 .
Figure 10.Gene signature risk score in the role of chemotherapy.(A) The correlation between different signature risk scores of OS patients and estimated IC50 values of bleomycin, cisplatin, doxorubicin, etoposide, and methotrexate in the TARGET (A), GSE21257 (B), and GSE39055 (C) datasets.The Wilcoxon signed-rank test was used to compare the estimated IC50 value of chemotherapy between groups.

Table 1 .
Characteristics of OS patients in the TARGET, GSE21257, GSE16091 and GSE39055 datasets.

Table 2 .
Prognostic genes obtained from the LASSO Cox regression model.

Table 3 .
Candidate research for comparison to our signature.

Table 4 .
Comparison of the nine-gene signature with previous published models in a univariate Cox analysis.CI confidence interval, HR hazard ratio, P P value of Cox regression model.The results were in bold if P < 0.05.